Factor analysis
library(readxl)
mydata <- read_xlsx("./Life values.xlsx")
mydata <- as.data.frame(mydata)
str(mydata)## 'data.frame': 18060 obs. of 18 variables:
## $ Cntry : chr "BG" "BG" "BG" "BG" ...
## $ Gndr : num 2 1 2 2 1 2 1 1 1 1 ...
## $ Age : num 76 43 50 51 70 31 40 48 71 71 ...
## $ X1_ipcrtiv : num 2 1 3 3 1 1 2 3 2 3 ...
## $ X2_ipeqopt : num 2 4 1 2 2 1 2 1 2 2 ...
## $ X3_ipshabt : num 5 2 4 3 1 1 1 1 1 2 ...
## $ X4_impdiff : num 1 1 3 3 2 3 3 3 1 4 ...
## $ X5_ipudrst : num 4 1 4 2 3 2 2 4 2 4 ...
## $ X6_ipmodst : num 4 6 4 6 6 1 2 3 2 3 ...
## $ X7_iphlppl : num 3 3 3 3 2 1 2 1 1 4 ...
## $ X8_ipsuces : num 5 1 4 4 2 1 2 2 2 5 ...
## $ X9_ipstrgv : num 3 2 6 1 1 1 1 1 2 2 ...
## $ X10_ipadvnt: num 3 1 2 4 1 5 3 3 4 6 ...
## $ X11_ipbhprp: num 2 5 4 2 5 1 1 2 2 3 ...
## $ X12_iprspot: num 5 1 4 4 5 3 2 2 2 3 ...
## $ X13_impenv : num 3 1 2 3 4 1 1 1 3 3 ...
## $ X14_imptrad: num 1 3 3 2 3 1 1 2 2 3 ...
## $ X15_impfun : num 5 1 4 4 1 3 4 2 3 6 ...table(mydata$Cntry)##
## BG CZ EE FI FR HR HU LT SI SK
## 2718 2476 1542 1577 1977 1592 1849 1659 1252 1418library(dplyr)##
## Attaching package: 'dplyr'## The following object is masked from 'package:psychTools':
##
## recode## The following object is masked from 'package:car':
##
## recode## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionmydata <- mydata %>%
filter(Cntry %in% c("SI")) #Selecting responses from SloveniaDescription
- cntry: Country
- gndr: Gender
- agea: Age of respondent, calculated
The following variables are measured on Likert scale:
1: Very much like me, 2: Like me, 3: Somewhat like me 4: A little like me, 5: Not like me, 6: Not like me at all
7: Refusal, 8: Don’t know, 9: No answer
- X1_ipcrtiv: Important to think new ideas and being creative
- X2_ipeqopt: Important that people are treated equally and have equal opportunities
- X3_ipshabt: Important to show abilities and be admired
- X4_impdiff: Important to try new and different things in life
- X5_ipudrst: Important to understand different people
- X6_ipmodst: Important to be humble and modest, not draw attention
- X7_iphlppl: Important to help people and care for others well-being
- X8_ipsuces: Important to be successful and that people recognise achievements
- X9_ipstrgv: Important that government is strong and ensures safety
- X10_ipadvnt: Important to seek adventures and have an exciting life
- X11_ipbhprp: Important to behave properly
- X12_iprspot: Important to get respect from others
- X13_impenv: Important to care for nature and environment
- X14_imptrad: Important to follow traditions and customs
- X15_impfun: Important to seek fun and things that give pleasure
# Identify the life-value indicators
value_items <- grep(
"^X[0-9]+_",
names(mydata),
value = TRUE
)
# Recode non-substantive responses only in these indicators
library(dplyr)
library(naniar)
mydata <- mydata %>%
mutate(
across(
all_of(value_items),
~ replace(.x, .x %in% c(7, 8, 9), NA)
)
)
# Retain complete responses on all life-value indicators
library(tidyr)##
## Attaching package: 'tidyr'## The following object is masked from 'package:reshape2':
##
## smithsmydata <- mydata %>%
drop_na(all_of(value_items))
# Create the dataset for factor analysis
mydata_FA <- mydata[value_items]FA
desc_stats <- data.frame(
N = apply(mydata_FA, 2, function(x) sum(!is.na(x))),
Missing = apply(mydata_FA, 2, function(x) sum(is.na(x))),
Mean = apply(mydata_FA, 2, mean, na.rm = TRUE),
Variance = apply(mydata_FA, 2, var, na.rm = TRUE),
Min = apply(mydata_FA, 2, min, na.rm = TRUE),
Max = apply(mydata_FA, 2, max, na.rm = TRUE)
)
round(desc_stats, 2)## N Missing Mean Variance Min Max
## X1_ipcrtiv 1178 0 2.21 1.14 1 6
## X2_ipeqopt 1178 0 1.70 0.62 1 6
## X3_ipshabt 1178 0 2.51 1.34 1 6
## X4_impdiff 1178 0 2.37 1.35 1 6
## X5_ipudrst 1178 0 2.04 0.79 1 6
## X6_ipmodst 1178 0 1.83 0.63 1 6
## X7_iphlppl 1178 0 1.81 0.56 1 5
## X8_ipsuces 1178 0 2.46 1.22 1 6
## X9_ipstrgv 1178 0 1.88 0.81 1 6
## X10_ipadvnt 1178 0 3.54 2.29 1 6
## X11_ipbhprp 1178 0 2.13 1.08 1 6
## X12_iprspot 1178 0 2.59 1.34 1 6
## X13_impenv 1178 0 1.72 0.56 1 6
## X14_imptrad 1178 0 2.42 1.33 1 6
## X15_impfun 1178 0 2.43 1.43 1 6R <- cor(mydata_FA) #Saving correlation matrixdet(R) #Determiant of correlation matrix## [1] 0.05797727library(psych)
cortest.bartlett(R, n = nrow(mydata_FA)) #Bartlett's test of sphericity## $chisq
## [1] 3335.136
##
## $p.value
## [1] 0
##
## $df
## [1] 105library(psych)
KMO(R) #KMO and MSA statistics## Kaiser-Meyer-Olkin factor adequacy
## Call: KMO(r = R)
## Overall MSA = 0.8
## MSA for each item =
## X1_ipcrtiv X2_ipeqopt X3_ipshabt X4_impdiff X5_ipudrst
## 0.84 0.78 0.77 0.81 0.82
## X6_ipmodst X7_iphlppl X8_ipsuces X9_ipstrgv X10_ipadvnt
## 0.81 0.86 0.79 0.84 0.76
## X11_ipbhprp X12_iprspot X13_impenv X14_imptrad X15_impfun
## 0.76 0.74 0.82 0.73 0.80library(psych)
fa.parallel(
mydata_FA,
fa = "fa",
fm = "minres",
cor = "cor",
n.iter = 1000,
sim = TRUE,
error.bars = TRUE
)
## Parallel analysis suggests that the number of factors = 4 and the number of components = NAlibrary(psych)
library(GPArotation)
factors <- fa_ordered(
mydata_FA,
nfactors = 4,
covar = FALSE,
fm = "minres",
rotate = "oblimin",
scores = "regression")
print.psych(
factors,
cut = 0.20,
sort = TRUE)## Factor Analysis using method = minres
## Call: psych::fa(r = r, nfactors = nfactors, rotate = "oblimin", scores = "regression",
## covar = FALSE, fm = "minres")
## Standardized loadings (pattern matrix) based upon correlation matrix
## item F1 F2 F3 F4 h2 u2 com
## X15_impfun 15 0.68 0.46 0.54 1.0
## X10_ipadvnt 10 0.61 0.43 0.57 1.1
## X4_impdiff 4 0.61 0.42 0.58 1.1
## X1_ipcrtiv 1 0.34 0.23 0.77 1.9
## X2_ipeqopt 2 0.60 0.36 0.64 1.1
## X5_ipudrst 5 0.54 0.33 0.67 1.1
## X13_impenv 13 0.41 0.22 0.28 0.72 1.7
## X7_iphlppl 7 0.41 0.25 0.35 0.65 2.0
## X6_ipmodst 6 0.36 0.32 0.29 0.71 2.2
## X14_imptrad 14 0.64 0.39 0.61 1.0
## X11_ipbhprp 11 0.52 0.38 0.62 1.4
## X9_ipstrgv 9 0.27 0.38 0.33 0.67 2.1
## X3_ipshabt 3 0.64 0.41 0.59 1.2
## X8_ipsuces 8 0.59 0.50 0.50 1.2
## X12_iprspot 12 0.27 0.54 0.40 0.60 1.6
##
## F1 F2 F3 F4
## SS loadings 1.55 1.45 1.29 1.24
## Proportion Var 0.10 0.10 0.09 0.08
## Cumulative Var 0.10 0.20 0.29 0.37
## Proportion Explained 0.28 0.26 0.23 0.22
## Cumulative Proportion 0.28 0.54 0.78 1.00
##
## With factor correlations of
## F1 F2 F3 F4
## F1 1.00 0.22 -0.02 0.45
## F2 0.22 1.00 0.30 0.19
## F3 -0.02 0.30 1.00 0.21
## F4 0.45 0.19 0.21 1.00
##
## Mean item complexity = 1.4
## Test of the hypothesis that 4 factors are sufficient.
##
## df null model = 105 with the objective function = 2.85 with Chi Square = 3335.14
## df of the model are 51 and the objective function was 0.17
##
## The root mean square of the residuals (RMSR) is 0.03
## The df corrected root mean square of the residuals is 0.04
##
## The harmonic n.obs is 1178 with the empirical chi square 79.08 with prob < 0.0071
## The total n.obs was 1178 with Likelihood Chi Square = 196.97 with prob < 5e-19
##
## Tucker Lewis Index of factoring reliability = 0.907
## RMSEA index = 0.049 and the 90 % confidence intervals are 0.042 0.057
## BIC = -163.68
## Fit based upon off diagonal values = 0.99
## Measures of factor score adequacy
## F1 F2 F3 F4
## Correlation of (regression) scores with factors 0.78 0.76 0.76 0.74
## Multiple R square of scores with factors 0.60 0.58 0.58 0.55
## Minimum correlation of possible factor scores 0.21 0.17 0.16 0.10library(dplyr)
mydata_FA <- mydata_FA %>%
select(-c("X1_ipcrtiv", "X6_ipmodst", "X9_ipstrgv"))# Recalculate the correlation matrix for the retained indicators
R <- cor(mydata_FA)
# Reassess sampling adequacy
psych::KMO(R)## Kaiser-Meyer-Olkin factor adequacy
## Call: psych::KMO(r = R)
## Overall MSA = 0.74
## MSA for each item =
## X2_ipeqopt X3_ipshabt X4_impdiff X5_ipudrst X7_iphlppl
## 0.72 0.74 0.80 0.78 0.80
## X8_ipsuces X10_ipadvnt X11_ipbhprp X12_iprspot X13_impenv
## 0.77 0.73 0.66 0.71 0.75
## X14_imptrad X15_impfun
## 0.63 0.77# Reassess the number of factors
psych::fa.parallel(
mydata_FA,
fa = "fa",
fm = "minres",
cor = "cor",
n.iter = 1000,
sim = TRUE,
error.bars = TRUE)
## Parallel analysis suggests that the number of factors = 4 and the number of components = NAlibrary(psych)
library(GPArotation)
factors <- fa_ordered(
mydata_FA,
nfactors = 4,
covar = FALSE,
fm = "minres",
rotate = "oblimin",
scores = "regression")
print.psych(
factors,
cut = 0.30,
sort = TRUE)## Factor Analysis using method = minres
## Call: psych::fa(r = r, nfactors = nfactors, rotate = "oblimin", scores = "regression",
## covar = FALSE, fm = "minres")
## Standardized loadings (pattern matrix) based upon correlation matrix
## item F1 F2 F3 F4 h2 u2 com
## X15_impfun 12 0.70 0.49 0.51 1.0
## X10_ipadvnt 7 0.64 0.47 0.53 1.1
## X4_impdiff 3 0.55 0.37 0.63 1.1
## X2_ipeqopt 1 0.60 0.36 0.64 1.1
## X5_ipudrst 4 0.54 0.31 0.69 1.1
## X13_impenv 10 0.46 0.30 0.70 1.5
## X7_iphlppl 5 0.43 0.32 0.68 1.7
## X3_ipshabt 2 0.62 0.41 0.59 1.2
## X8_ipsuces 6 0.58 0.48 0.52 1.2
## X12_iprspot 9 0.58 0.42 0.58 1.4
## X14_imptrad 11 0.66 0.43 0.57 1.0
## X11_ipbhprp 8 0.46 0.34 0.66 1.8
##
## F1 F2 F3 F4
## SS loadings 1.38 1.22 1.19 0.91
## Proportion Var 0.12 0.10 0.10 0.08
## Cumulative Var 0.12 0.22 0.32 0.39
## Proportion Explained 0.29 0.26 0.25 0.19
## Cumulative Proportion 0.29 0.55 0.81 1.00
##
## With factor correlations of
## F1 F2 F3 F4
## F1 1.00 0.25 0.44 -0.03
## F2 0.25 1.00 0.21 0.22
## F3 0.44 0.21 1.00 0.20
## F4 -0.03 0.22 0.20 1.00
##
## Mean item complexity = 1.3
## Test of the hypothesis that 4 factors are sufficient.
##
## df null model = 66 with the objective function = 2.09 with Chi Square = 2446.14
## df of the model are 24 and the objective function was 0.1
##
## The root mean square of the residuals (RMSR) is 0.02
## The df corrected root mean square of the residuals is 0.04
##
## The harmonic n.obs is 1178 with the empirical chi square 42.22 with prob < 0.012
## The total n.obs was 1178 with Likelihood Chi Square = 115.59 with prob < 5.9e-14
##
## Tucker Lewis Index of factoring reliability = 0.894
## RMSEA index = 0.057 and the 90 % confidence intervals are 0.047 0.068
## BIC = -54.13
## Fit based upon off diagonal values = 0.99
## Measures of factor score adequacy
## F1 F2 F3 F4
## Correlation of (regression) scores with factors 0.77 0.76 0.74 0.74
## Multiple R square of scores with factors 0.60 0.57 0.55 0.55
## Minimum correlation of possible factor scores 0.19 0.14 0.10 0.10X15_impfun: Important to seek fun and things that give pleasure
X10_ipadvnt: Important to seek adventures and have an exciting life
X4_impdiff: Important to try new and different things in life
X1_ipcrtiv: Important to think new ideas and being creative
X2_ipeqopt: Important that people are treated equally and have equal opportunities
X5_ipudrst: Important to understand different people
X13_impenv: Important to care for nature and environment
X7_iphlppl: Important to help people and care for others well-being
X3_ipshabt: Important to show abilities and be admired
X8_ipsuces: Important to be successful and that people recognise achievements
X12_iprspot: Important to get respect from others
X14_imptrad: Important to follow traditions and customs
X11_ipbhprp: Important to behave properly
indicator_labels <- c(
X1_ipcrtiv = "Important to think new ideas and being creative",
X2_ipeqopt = "Important that people are treated equally and have equal opportunities",
X3_ipshabt = "Important to show abilities and be admired",
X4_impdiff = "Important to try new and different things in life",
X5_ipudrst = "Important to understand different people",
X6_ipmodst = "Important to be humble and modest, not draw attention",
X7_iphlppl = "Important to help people and care for others well-being",
X8_ipsuces = "Important to be successful and that people recognise achievements",
X9_ipstrgv = "Important that government is strong and ensures safety",
X10_ipadvnt = "Important to seek adventures and have an exciting life",
X11_ipbhprp = "Important to behave properly",
X12_iprspot = "Important to get respect from others",
X13_impenv = "Important to care for nature and environment",
X14_imptrad = "Important to follow traditions and customs",
X15_impfun = "Important to seek fun and things that give pleasure")
factor_labels <- c(
F1 = "Stimulation",
F2 = "Self-transcendence",
F3 = "Achievement",
F4 = "Tradition and Conformity")
loadings_diagram <- as.matrix(factors$loadings)
rownames(loadings_diagram) <- indicator_labels[rownames(loadings_diagram)]
colnames(loadings_diagram) <- factor_labels[colnames(loadings_diagram)]
Phi_diagram <- factors$Phi
dimnames(Phi_diagram) <- list(colnames(loadings_diagram), colnames(loadings_diagram))
op <- par(cex = 1.5)
psych::fa.diagram(
loadings_diagram,
Phi = Phi_diagram,
cut = 0.40,
simple = TRUE,
main = "Factor Analysis")
par(op) Residuals
residuals <- factors$residual #Residuals
residuals <- as.matrix(residuals[upper.tri(residuals)]) #Only taking into account one half of residual matrix
highresiduals <- abs(residuals) > 0.05 #Logical test for finding high residuals
sum(highresiduals)/nrow(residuals) #Share of high residuals## [1] 0.04545455Saving factor values
head(factors$scores) #Showing factor values## F1 F2 F3 F4
## [1,] -0.4002379 -1.257436721 -0.8889618 -0.5139052
## [2,] 0.6692058 0.511099087 1.2064647 0.2747212
## [3,] 1.4259899 -0.689320740 0.3276818 -0.4914667
## [4,] -0.5490813 -0.001598386 -0.3502387 0.2328948
## [5,] 1.1822772 -1.247398156 -0.8811219 0.4444021
## [6,] 0.7697104 -0.863318982 -0.6932335 -1.1231987scores_ordered <- factors$scores
mydata$F1_Stimulation <- scores_ordered[, "F1"]
mydata$F2_Self_transcendence <- scores_ordered[, "F2"]
mydata$F3_Achievement <- scores_ordered[, "F3"]
mydata$F4_Tradition_and_Conformity <- scores_ordered[, "F4"]